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The magnetization switching dynamics of biaxial single-domain homogeneous ferromagnetic parti- 
cles, in which the two easy axes are perpendicular to each other, is simulated using a 4-state clock 
model. A zero-field mapping of the statics between the symmetric 4-state clock model and two 
decoupled Ising models is extended to non-zero field statics and to the dynamics. This significantly 
simplifies the analysis of the simulation results. We measure the magnetization switching time of 
the model and analyze the results using droplet theory. The switching dynamics in the asymmetric 
model is more complicated. If the easy axis is perpendicular to the stable magnetization direction, 
the system can switch its magnetization via two different channels, one very fast and the other 
very slow. A maximum value for the switching field as a function of system size is obtained. The 
asymmetry affects the switching fields differently, depending on whether the switching involves one 
single droplet or many droplets of spins in the stable magnetization configuration. The angular 
dependence of the switching field in symmetric and asymmetric models is also studied. 

PACS Number(s): 75.40.Mg, 05.50.+q, 02.70.Lq, 02.50.-r 
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I. INTRODUCTION 

Biaxial media have recently attracted much research 
interest as potential materials for applications in high- 
density magnetic recording. This is a consequence of 
their good thermal stability, high coercivity and coer- 
civity squareness, and low medium noise (!]]- Q. With 
the current rapid increase in recording density, the size 
of a magnetic particle used to store a single bit may soon 
become so small that it can only contain one single mag- 
netic domain. Understanding the dynamics of magneti- 
zation switching in individual single-domain particles is 
therefore essential. Biaxial media may be obtained, for 
example, by sputtering a Co-based alloy onto a Cr sub- 
strate by fabricating perovskite superlattices con- 
sisting of different ferromagnetic metallic oxides 0, or 
by epitaxially growing Fe on GaAs(OOl) ||. By apply- 
ing lithographic patterning techniques to precisely con- 
trol the particle shape, a uniaxial medium may also ex- 
hibit biaxial properties due to competition between the 
magneto-crystalline and shape anisotropics pJlOJ. The 
two easy axes in biaxial media are often perpendicular to 
each other, and they may have the same or different mag- 
netic anisotropies. In thin films, these two axes can both 
be in-plane, or one of them can be out-of-plane, often 
perpendicular to the plane of the film. New techniques, 
such as nanolithography and ultra-high resolution scan- 
ning microscopies, have recently made the synthesis and 
experimental observation of isolated single-domain parti- 
cles possible. For instance, the magnetization switching 
of individual barium ferrite nanoscale particles has been 
observed using Magnetic Force Microscopy (MFM), and 
a maximum switching field was observed for particles of 



diameter near 55 nm For other particles, the an- 

gular dependence of the switching (or coercive) field has 
been measured [[l2 13 , the switching dynamics investi- 
gated 12 1, and the switching statistics measured [fl3|| . 
In Refs. [14|] and (l5), arrays of single-domain Co dots 
were fabricated, and the magnetization switching was ob- 
served by MFM. Depending on the shape of the Co dot, 
both in-plane and out-of-plane magnetic moments were 
observed. 

A powerful non-perturbative theoretical method to 
study models of magnetization switching dynamics is 
Monte Carlo simulation. One standard model used to 
study uniaxial magnets is the kinetic Ising model [jl6],|l7j • 
However, the Ising model can only be used as a reason- 
ably realistic representation of extremely anisotropic uni- 
axial magnets, and a different model is required to cap- 
ture the unique properties of biaxial materials. In this 
work we study a 4-state clock model with two perpen- 
dicular axes. In the most symmetric case there is no 
applied field and four states are equivalent. Suzuki has 
found a mapping between this symmetric 4-state clock 
model and two superimposed Ising models |Q. When a 
field is applied, this symmetry is lost and o the results of 
extensive simulations of the switching dynamics. In or- 
der to simplify the analysis of the simulation results, we 
extend Suzuki's mapping to non-zero applied field and to 
the dynamics. 

The rest of this paper is organized as follows. In Sec 
we describe the model used in the simulations. In Sec. 
we briefly review the aspects of droplet theory used 
analyze the numerical results. In Sec. IV and Sec. 
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we present and analyze our numerical results for lifetime 
and switching field, respectively. Finally, conclusions are 
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drawn in Sec. V] 



II. MODEL DESCRIPTION 

The 4-state clock model is denned by the microscopic 
Hamiltonian 

H = — J pj- pj — H ■ (2j — Ay]{pi ■ e ) 2 ■ (1) 

(id) i 1 

The unit vector pi is the orientation of the spin at site 
i, which can be in one of 4 directions: eb, ex, e*2, and 
e*3, as shown in Fig. [jj The angle between the applied 
magnetic field H and eb is 9. The dimensionless system 
magnetization is given by 



m = L 



(2) 



We consider a two-dimensional square lattice with L 
spins on each side, and periodic boundary conditions are 
applied in both directions to suppress boundary effects. 
(Boundary effects for the fsing model are discussed in 
Rcf. pS% .) The lattice constant is set to unity. In Eq. (|l]), 
J2i runs over all L? sites, and Ylu j) runs over au nearest- 
neighbor pairs on the lattice. The first term in Eq. (Q) is 
the energy due to exchange interactions between nearest- 
neighbor spins, and J > is the ferromagnetic exchange 
interaction. The second term is the energy, L 2 H ■ rh due 
to the applied field H . The third term represents the 
energy due to a different preference for one of the two 
axes, eb-eb and ei-eV If A = 0, the model has 4-fold 
symmetry for H — 0, and the statics of this model is 
equivalent to two decoupled Ising models Jl8| |. In Ap- 
pendix A, this mapping is extended to non-zero fields 
and to the dynamics. A nonzero value of A destroys this 
4-fold symmetry. Therefore, throughout this paper we 
refer to the model with A ^ as asymmetric. (In the 
statistical mechanics literature, clock models with A ^ 
are often referred to as "anisotropic." However, in this 
paper we use the terms "anisotropy" and "anisotropic" 
as they are customarily understood in the literature on 
magnetic materials.) For A > and H = 0, the eb- 
e*2 axis is favorable energetically, while the ei-eb axis is 
preferred for A < 0. 




The critical temperature, T c , of the model with A = 
and H — is known exactly. Its value is k-gT c /J 
1 1m; 1 + y/2) w 1.135 @, exactly that of a two- 
dimensional square-lattice Ising model with interaction 
constant J/2. Here fee is Boltzmann's constant. The 
simulation temperature is chosen to be T — 0.8T C . This 
simulation temperature is high enough to get a reason- 
able Monte Carlo acceptance rate, while it is low enough 
to avoid complications near the critical point such as crit- 
ical slowing-down and finite-size effects J2Q] . Based on 
studies of the kinetic Ising model pq , pi|,p2[ , simulation 
results at other moderately low temperatures should be 
qualitatively similar to those obtained here. 

The initial configuration of the system is ordered, with 
pi = eo for all spins. This may be achieved by apply- 
ing a strong field in the eb direction for time t < 0. 
Then at t = the external magnetic field is changed 
to H (Fig. Q). During the simulation, first one spin is 
selected at random, then a new candidate state is cho- 
sen. At low temperatures, large orientational changes of 
the spin (eb <-* b\ and e*i <-> eb) are very unlikely [p3| . 
Consequently, the candidate state is obtained by rotat- 
ing the spin by 90° or —90°. The a priori probabilities 
of clockwise and anti-clockwise rotations are each set to 
1/2. 

The probability that the spin flips to the chosen state 
is given by Glauber dynamics |24| 



W 



exp(-/3A£:) 
1 + cxp(-/?A£0 ' 



(3) 



where AE gives the change in the energy of the system 
that would result if the spin flip were accepted, and /3 = 
l/k^T. The Glauber dynamic is chosen because it can be 
derived from quantum mechanics under certain stringent 
conditions [^5| . Here time is reported in the unit of one 
Monte Carlo steps per spin (MCSS). 

The lifetime is defined as the average time, (r), needed 
to satisfy a specific stopping criterion (discussed below), 
where r is the first-passage time to the stopping criterion 
for an individual Monte Carlo run. The average is taken 
over at least 1000 independent Monte Carlo runs in order 
to obtain good statistics. In addition to the mean, we also 
measure the relative standard deviation 



r = ^(r 2 ) - (r)V(r) 



(4) 



FIG. 1. Vector representation of the 4-state clock model. 



The switching field is measured by observing whether 
on average the system magnetization has reversed after 
a fixed waiting time, here chosen as t w = 200 MCSS. 
As Ising model studies have demonstrated pq , pl| , the 
switching fields for different waiting times are qualita- 
tively similar. 

To measure the lifetime of a metastable state, a crite- 
rion for escape from the metastable state must be chosen. 
In studies of the kinetic Ising model, the usual criterion is 
that an escape is registered and the simulation stopped 
the first time the magnetization reaches a fixed cutoff 



2 



value, m s to P - (In other words: the stopping time is the 
first-passage time to m s t op -) Results for different values 
of m stop are significantly different only in the determin- 
istic regime, and even in this regime the results are qual- 
itatively similar p6| . 

In some experimental studies of magnetization switch- 
ing, the sign of the magnetization along a particular 
direction is much easier to measure than its magni- 
tude pi 12\. The stopping criterion most directly rel- 
evant to these experiments is therefore to use m s t p = 0. 
For the Ising model, this is obtained by stopping the 
simulation when half of the spins have entered the sta- 
ble state (and half have left the metastable state). For 
other experimental situations, a different criterion might 
be easier to measure. However the stopping criterion 
chosen should not have a large effect on the results. 

In order to completely map the the 4-state clock model 
onto two superimposed Ising systems (see Appendix A), 
the stopping criterion for the 4-state clock model must 
map onto the stopping criteria for the two individual 
Ising systems. If the stopping criteria for both Ising 
systems are given by rn st op = 0, two different stopping 
criteria for the 4-state clock model can be naturally en- 
visioned. Denoting these by S u and Sn, they are given 
in terms of the first passage to the stopping parameter 
Mstop by 



Su has 
Sn has/Li; 



stop 



Mstop = m I0 m n 
- 1 - ry(-m/o)?y(-?Ti/i) 



(5) 
(G) 



where each toj q for a = 0, 1 is the time-dependent mag- 
netization of one of the Ising systems and rj(x) is the 
unit step function. The criterion Su is to stop when the 
first Ising system reaches zero magnetization, while for 
Sn the simulation is continued until the last of the two 
Ising systems reaches zero magnetization. 

The mapping between the clock model and the two 
Ising models illustrates that in order for a clock spin to 
escape from the initial state, only one of the two Ising 
spins on the site needs to flip, while both Ising spins on 
the same site must flip in order for the clock spin to enter 
the stable state. Therefore, for the clock model the Su 
stopping criterion implies that > ^ of the clock spins have 
escaped from the initial state. The Sn stopping criterion 
means a certain fraction (< 5) of the clock spins have 
entered the stable state. For the symmetric clock model, 
the two Ising models are decoupled, and their switching 
processes are independent of each other. Therefore, for 
S n on average 1/4 of the clock spins will have entered 
the stable state at the time the simulation is stopped. In 
this paper, we use only the Sn stopping criterion since it 
more closely models an experimental stopping criterion 
of zero magnetization. For the asymmetric model, we 
also use this stopping criterion in order to observe the 
effects of A and H consistently. In this case the average 
fraction of clock spins in the stable state at the end of 
the simulation is between 1/4 and 1/2, depending on A. 
A more detailed discussion of the quantitative differences 



between different stopping criteria is given in Rcf. 



III. DROPLET THEORY 

During magnetization switching in single-domain par- 
ticles, small "droplets" of the stable phase are continually 
created and destroyed by thermal fluctuations [^6[ . The 
fate of a droplet is determined by the competition be- 
tween its total surface free energy, which is proportional 
to its surface area, and the coupling of its magnetiza- 
tion with the applied field, which is proportional to its 
volume. If the droplet is small, further growth produces 
an increase of the free energy due to the surface tension. 
Therefore, the droplet will be destroyed with high proba- 
bility. If the droplet is sufficiently large, this free-energy 
penalty is offset by the benefit obtained from increas- 
ing the number of spins parallel to the applied magnetic 
field. Therefore the droplet is very likely to grow further. 
For a droplet with a critical radius R c , the tendencies 
towards growth and shrinkage balance, and the droplet 
will grow or shrink with equal probability. By comparing 
R c with the system size L, the lattice spacing a, and the 
average distance between supercritical droplets Rq, the 
magnetization switching dynamics is divided into four 
regimes ]l6|,|l]]26|||] : 

• The Coexistence (CE) regime [a < L < R c ], in 
which the switching is governed by subcritical fluc- 
tuations on the scale of the system size. 



• The Single-droplet (SD) regime [a < R c < L < 
i?o] ! characterized by switching via a single critical 
droplet. Once a supercritical droplet is nucleated, 
it will most likely grow to occupy the entire system 
before another droplet is nucleated. 

• The Multi-droplet (MD) regime [a < R c < R < 
L], characterized by switching via a finite density 
of critical droplets. 

• The Strong-field (SF) regime [a ~ R c ~ Rq < L], 
in which the concept of a droplet is no longer ap- 
plicable. 

Since the symmetric 4-state clock model is equivalent 
to two decoupled Ising models, we can simplify the anal- 
ysis of the clock model based on extensive studies of ki- 
netic Ising models |l^ , ^l| , ^8| , p6| . Here a brief review of 
the droplet theory for the Ising model is given. For more 
details, see the references listed above. 

In the SD regime, the lifetime is given by 



n » [L d r(T, Hi) 



(7) 



where Hi is the field in the Ising Hamiltonian, and the 
nucleation rate per unit volume T(T, H\) is pl| 



T{T,Hi)*iA(T)\Hi\ K ex9 



0(H?)) 



(8) 
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Here A(T) is a non-universal prefactor, and K is believed 
to be 3 for two dimensional (d = 2) Ising systems |2l]] . At 
T = 0.8T C , 3(T)/ J? « 0.92 @], where Ji is the exchange 
interaction in the Ising Hamiltonian. 

In the MD regime, the system can be approximately 
partitioned into (L/Ro) d ^> 1 cells, each of which de- 
cays via an independent Poisson process with rate RqT. 
Therefore, the volume fraction is self-averaging, and the 
average lifetime (n), defined as the average first-passage 
time to zero magnetization, is given by the well-known 
Avrami theory of metastable decay as p6|.p6| 



(ti) a Offerer, Bk)] 



-i/(d+i) 



(9) 



By taking the derivative of the logarithm of (ti) with 
respect to | Ji/Hi]^ 1 , one obtains 



A 



icir 



dln(n) 



" d(\Jl/Hl\ d ^) 



A|ffi/Ji 



id— 1 



A 



J 



d-l 



(10) 



with A = K/(d - 1) and A = (3E(T) in the SD regime, 

while A = (K + d)lid 2 ~ 1) and A = E i T )/i d + 1) in 
the MD regime 

For the symmetric 4-state clock model, if the applied 
field is along the e*2 direction, the applied field H\ and 
the exchange interaction Ji of the Ising model are related 
to the corresponding parameters for the clock model by 
Hi = H/2 and Ji = J/2 (Table [|. By replacing the 
parameters for the Ising model by the corresponding ones 
for the clock model in the above analysis, we get the 
dimensionless value A e g/ J d_1 for the clock model, which 
is 



A 



cir 



dln(- 



A 



3 d - 1 d{\J/H\ d - 1 ) 



I^=X\H/J\ d ^ + IJW ^, (11) 



Since Ji = 1/2 J, comparing Eq. ( |i~0| ) and Eq. (pi]), we 

have 



A ff / J d 1 — Ai C ff / J d 1 



(12) 



Therefore we can compare the simulation results for the 
symmetric 4-state clock model directly with the theoret- 
ical predictions for the Ising models. For the asymmetric 
clock model, we calculate the same derivative. 



implies that both constituent Ising lattices must have 
reached zero magnetization. Therefore, the lifetime of 
the clock model is the longer of the lifetimes of the two 
Ising models. If the applied field is in the e*2 direction, 
the two Ising systems are identical, and the probability 
density of the lifetime for the clock model can be simply 
written as 

Pc (t) = (n)- 1 exp(-r/2(7i))[l - exp(-r/2(n))] . (14) 

The mean and relative standard deviation of this distri- 
bution are (t) = 3(ti) and r = \/5/3, respectively. The 
derivation of p c (r) for the general case that H is not par- 
allel to e*2 is detailed in Appendix B with the general 
result given as Eq. (pfj). 

In the MD regime, since the switching process of the 
Ising model is deterministic, based on the mapping be- 
tween the Ising model and the 4-state clock model, the 
switching dynamics of the 4-state clock model is also de- 
terministic. By "deterministic" we mean that the life- 
time distribution is characterized by a very small rela- 
tive standard deviation, r <C 1. Figure [2] and Fig. || 
show (r) and r as functions of the inverse applied field, 
J/\H\, for different lattice sizes. In these measurements, 
the applied field H is always parallel to €2 ■ In the panels 
showing r, the horizontal solid and dashed lines corre- 
spond to r = 1 and r = v5/3, respectively. As pre- 
dicted, for the isotropic clock model, r approaches \/5/3 
as \H\ decreases. Also shown in these figures are results 
for the asymmetric clock model. When the easy axis 
is along e"o— §2 (A > 0), the lifetime is longer than that 
of the symmetric model with the same applied field be- 
cause a stronger field is needed to flip the spins from the 
metastable state eb to the intermediate states, e\ and e*3. 
When the easy axis is along e\-e^ (A < 0), however, the 
lifetime is shorter than that for the symmetric model in 
the MD regime, and longer in the SD regime. The life- 
time distribution is wider in the SD regime, and r has 
a system size-dependent maximum for A < 0. In either 
CclSG, clS I if I decreases, the asymmetry dominates, and 
the two states along the easy axis are strongly favored 
over the other two states. Therefore for weak H we ex- 
pect the asymmetric clock model to behave essentially 
like a single Ising model, so that r will approach unity. 



IV. LIFETIME MEASUREMENTS 



The statistical properties of the lifetime depend on the 
switching regime of the system. For an Ising system in 
the SD regime, the probability density of the lifetime is 
exponential ]T^ , ^l| , p6| , p8| 

p I (r) = (r I )" 1 exp(-r/(r I )), (13) 

where (ti) is the mean lifetime. As discussed in Sec. |J, 
for the 4-state clock model, magnetization switching, as 
defined by the stopping criterion used in this paper, 
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FIG. 2. (r) and r vs. J/\H\ for L = 100. 
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FIG. 3. (r) and r «s. J/|fl"| for L = 32. 

In order to understand the switching process for the 
asymmetric model better, we show selected snapshots 
of the model during the simulations in Figs. ^-|[ The 
darkest regions represent the initial state, the lightest 
ones represent the stable state, and the two intermedi- 
ate shades of gray represent the intermediate states. The 
applied field is chosen so that r is near its maximum for 
A < (Fig.|). For A > (Fig. §), there is only one 
channel for the model to switch the magnetization: first 
a droplet of one intermediate state is nucleated and grows 
quickly. Then a droplet of the stable state nucleates and 
grows on this intermediate background. This process is 
similar to single-droplet decay in the Ising model, except 
that there is no intermediate phase in the Ising model. 
For A < 0, however, two channels exist. In the slow 
channel (Fig. ^|), the system first enters a uniform inter- 
mediate state, from which it takes a long time to nucleate 
a stable droplet. In the fast channel (Fig. |J), droplets in 
different intermediate states are formed simultaneously. 
As they meet, a stable droplet is quickly nucleated on 
their common boundary to decrease the free energy. Be- 
cause of the existence of these two channels, the relative 
standard deviation of the lifetime becomes very large, ex- 
ceeding unity, which is the maximum value for a single 
Ising model. 
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FIG. 4. One typical switching event for A/ J = 0.025, 
L = 100, \H\/J = 0.09. r = 66000 MCSS. 
(r) = 25886 MCSS, r = 0.82. The times for these four snap 
shots are t = 62010, 62681, 65130, and 65980 MCSS, respec- 
tively. 
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FIG. 6. One typical realization of decay via the fast channel 
for A/ J = -0.025, L = 100, \H\/J = 0.09. r = 1649 MCSS. 
(r) = 5659 MCSS, r = 3.02. The times for these four snap 
shots are t = 410, 820, 1230, and 1640 MCSS, respectively. 

Figure [7] shows the derivatives of the logarithm of the 
lifetime with respect to J/\H\, A e g of Eq. (|l0|), for differ- 
ent system sizes. In this figure, the straight lines are the 
theoretical results obtained from Eq. (|o|). Indeed, for 
the symmetric clock model, our numerical results show 
good agreement with the theory. Note that there are no 
adjustable parameters used in obtaining the theoretical 
results. In the asymmetric cases, when the applied field 
is strong the effect of the asymmetry is small. For weak 
fields, the asymmetry dominates the switching process, 
and the results for the clock model cannot be mapped 
onto those for the Ising model. Unfortunately, we could 
not extend the simulations deep into this regime, due to 
the prohibitive computer resources required: it takes at 
least two months to simulate 1000 individual realizations 
for L = 100 and J/\H\ > 12.5 (Fig. |) on a single node 
of an IBM sp2 computer. 



FIG. 5. One typical realization of decay via the slow 
channel for A/J = -0.025, L = 100, \H\/J = 0.09. 
r = 44542 MCSS. (r) = 5659 MCSS, r = 3.02. The times 
for these four snap shots are t = 1100, 39960, 43130, and 
44400 MCSS, respectively. 
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FIG. 7. Aeff/J fs. |H|/J. The solid lines are the theoreti- 
cal predictions for the Ising model from Eq. (|ic|). 

V. SWITCHING-FIELD MEASUREMENTS 

We define the switching field as follows. We apply 
a field H to the system in the initial state for a cho- 
sen waiting time t w = 200 MCSS. If the magnetization 
has switched N/2 times in N simulations, then the mag- 
nitude of the applied field is called the switching field, 
H svr (t w ). We choose this definition in order to simulate 
qualitatively the experiment discussed in Ref. 0|. The 
relationship between the switching field and the coerciv- 
ity is that in coercivity measurements, the applied field 
is changed sinusoidally or linearly p9| , while in our sim- 
ulation the applied field is reversed abruptly. 

First we compare the switching field for the symmetric 
clock model and the Ising model. The switching field as 
a function of the system size L is illustrated in Fig. [|. 
Since the time scale of the clock model is twice that for 
the Ising model, the switching field for the Ising model 
is measured with waiting time t w — 100 MCSS. In the 
MD regime, the lifetime for the symmetric clock model is 
twice that for the Ising model, thus the switching fields 



for the clock model and the Ising model are the same. 
In the stochastic regime, from Eq. (|l|) and Eq. (Q, the 
lifetime for the clock model is three times that of the 
Ising model. Therefore, the switching field for the clock 
model is larger than that for the Ising model. Here we 
emphasize that the maximum shown in Fig. |^ is the result 
of different switching mechanisms dominating the decay 
in the CE, SD and MD regimes. It is not the effect of 
demagnetization since there are no dipole-dipole interac- 
tions in our Hamiltonian. For more detailed discussions 
of the dependence of H sw on L, t w and T in the Ising 
model, see Ref. Q. 
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FIG. 8. Switching fields for the symmetric clock model and 
the Ising model. 

Next we consider the effects of the asymmetry A on the 
switching field. The results are shown in Fig. ^. Again, 
H is in the direction. When the easy axis is along e$- 
e*2 (A > 0), increasing A leads to a monotonic increase 
of the switching field. When A < 0, on the other hand, 
making A more negative decreases the switching field in 
the MD regime. However, in the SD regime different A 
shows little effect on the switching field. 
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FIG. 9. Effect of A on H BW for 9 = 180°. 
The effect of the direction of the applied field on the 
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switching field is shown in Fig. [lO|. Here H sw is shown as 
a function of L for different 9 at A = 0. The symmetric 
clock model is equivalent to two decoupled Ising models, 
and the field components acting on these two Ising sys- 
tems change with 8. Thus the switching behaviors of the 
two Ising systems are different, and the switching of the 
clock model is determined by the slower one. Therefore 
H sw increases as H deviates further away from the e~2 
direction. 
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FIG. 11. 0- and A- dependence of the switching field for 
L = 8 and 40. 



VI. CONCLUSION 



FIG. 10. Effect of on # sw for A = 0. 

Finally we have measured the ^-dependence of H sw for 
two system sizes (L = 8 and L = 40) with different val- 
ues of A. The results are shown in Fig. O. When the 
direction of H is close to e*2 , the energetic preference for 
the two intermediate states is similar and mainly deter- 
mined by the asymmetry. Therefore H sw is governed by 
A. As 9 approaches 135°, H increases the stability of one 
intermediate state and decreases that of the stable state. 
When A < 0, this tendency is further amplified, there- 
fore the switching field increases. For A > 0, H makes 
the flip from eo to the favored intermediate state easier, 
whereafter A acts like an applied field in flipping the spin 
to e*2. Thus H sw decreases compared with the isotropic 
model. 
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In this work we have studied magnetization switch- 
ing in biaxial magnetic media by large-scale Monte Carlo 
simulation and theoretical analysis of a kinetic 4-state 
clock model. The model reduces the infinite anisotropy 
in the kinetic Ising model by introducing two intermedi- 
ate magnetization directions between the magnetization 
states parallel and antiparallcl to the applied field. If the 
two perpendicular axes are equivalent, the switching pro- 
cess is equivalent to that of two decoupled Ising models. 
With asymmetry, the clock model shows a more com- 
plicated switching behavior. When the axis along the 
initial state is favored (A > 0), the model becomes Ising- 
like. On the other hand, if the other easy axis is favored, 
the system can switch its magnetization via two different 
channels, which can lead to a wider distribution of the 
lifetimes in the SD regime. Another factor affecting the 
switching process is the direction of the applied field. We 
studied the direction dependence of the switching field in 
the isotropic model. With the applied field deviating 
further from the stable state, a stronger switching field is 
needed. The different switching mechanisms in the CE, 
SD, and MD regimes lead to a maximum switching field 
for a specific system size, in agreement with experimen- 
tal results. We conclude that the 4-state clock model is 
a simple but informative platform to study the magneti- 
zation switching in biaxial magnetic particles. 
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APPENDIX A 

In this appendix, the dynamical mapping between the 
Ising model and the symmetric 4-state clock model is 
derived. 

First we define two unit vectors, so and si, such that 
the angle between s (si) and e is -45° (+45°) (Fig.[l2|). 
Let fli denote the direction of the clock spin at site i, 
which can be eo, ei, €2, or S3. The unit vector /7j can be 
expressed in terms of Sq and si as 



Mi 



1 
75 



(15) 



where ujoi and uju can take the values ±1. The applied 
field can also be decomposed as 



where 



H = Hos + HiSi 



H = \H\ cos(45° + 6) 
Hi = \H\ sin(45° + 6) 



(16) 

(17) 
(18) 



sum of two decoupled Ising Hamiltonians with exchange 
interactions J/2. The applied field in each Ising Hamilto- 
nian equals the field component times l/y/2. Therefore, 
the 4-state clock model in equilibrium can be treated as 
two super imp osed Ising systems, with two Ising spins at 
each site JlSj . 

The next step is to construct a mapping of the spin-flip 
dynamics between the 4-state clock model and the Ising 
models. Since each spin can only flip to its neighboring 
state, this can be done by choosing one of the two Ising 
systems randomly and flipping the corresponding Ising 
spin. Because each time only one Ising lattice is selected, 
the change of the energy due to the spin flip is obtained 
by updating either Ho or Hi in Eq. (fL9[). 

With the above analysis, it follows that the dynamics 
of the symmetric 4-state clock model is equivalent to that 
of two decoupled Ising models. Not only are the spin- 
flip mechanisms interchangeable, but also only one Ising 
Hamiltonian in Eq. ( |l9|) is affected in each flip because 
Ho and Hi are decoupled. Since each clock spin can 
be thought to be composed of two Ising spins, during 
the simulation the time scale for the 4-state clock model 
is twice that for the Ising model. The mapping of the 
parameters between the symmetric 4-state clock model 
and Ising models is summarized in Table B In this table, 
t is the time measured in Monte Carlo steps per spin 
(MCSS). 

TABLE I. Mapping between the symmetric 4-state Clock 
and Ising Models 



Clock 


Ising system #0 


Ising system #1 


J 


J/2 


J/2 


H 


^|tf|cos(45 +9) 


i|i?|si n (45°+0) 


t 


t/2 


t/2 



ex. 



H *, 



FIG. 12. Mapping of the 4-state clock model onto a super- 
position of two Ising models. 



Substituting Eqs. (|T|)-@ into Eq. (Q), we get 



H = Ho + Hi 



LOOi^Oj 



UJliLOij 



Ho \ - 

Hi ^ 



(19) 
(20) 

(21) 



From these equations, it is clear that the Hamiltonian of 
the symmetric 4-state clock model is equivalent to the 
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APPENDIX B 

In this appendix, the lifetime distribution for the 
isotropic 4-state clock model is derived. We denote the 
lifetimes of the two Ising models by To and ti, and the 
lifetime of the 4-state clock model by r c . With the stop- 
ping criterion used in this work, we have 



Pc(t) 



t c = 2max(T ,Ti) 



(22) 



Thus, for a given value r of the lifetime, the probability 
that t c < t can be written as 



P(t c <t) = P(max(ro,n) < t/2) 



(23) 



Here the factor 2 comes from the mapping of the time 
scale between the clock and the Ising models given in 
Table D. Since the two Ising models are decoupled, the 
lifetime distributions are independent. Therefore 



P(r c < r) = P(r < t/2) ■ P( n < t/2) 



(24) 



In the SD regime, the lifetime of each Ising model is 
exponentially distributed, 



(25) 



P(r k < t) = 1 - exp(-T/(r fc )) , 



where k = 0, 1. Substituting Eq. (]25|) into Eq. (g4j) and 
differentiating with respect to r, we obtain the probabil- 
ity density of the lifetime of the clock model as 



Pc(t) 



exp(-T/2(r )) [1 - exp(-r/2(n))] 
2<ro> 

exp(-r/2(r 1 )) [1 - exp(-r/2(r ))] 

2(n) 



(26) 



If (t > = (n) = fa), we get Eq. <y>. 

In the MD regime, since the magnetic switching is de- 
terministic, the lifetime probability density of the Ising 
model is well approximated by a (5-function. Thus for the 
Ising model, 



P(ri < r) s» t){t - ti) , 



(27) 



where rj is the unit step function. Following the same 
procedure as above we get the lifetime distribution for 
the clock model, 



P c (t) « 6(t - 2max[(r ), (n)]) 



(28) 



Finally we consider the case that the field is applied 
in such a direction that one Ising system is in the MD 
regime with average lifetime (tmd), and the other in the 
SD region with the average lifetime (tsd)- Starting from 
Eq. (p4|), we then have 

Pc(r c <t)= r?(r - 2r MD )[l - exp(-r/2(r SD ))] , (29) 

which yields 



[1 - exp(-(TMD)/2(r SD ))] ( 5(T - 2(r MD )) 
exp(-T/2(r SD )) 



2(t sd ) 



-t}(t - 2(t M d)) 



(30) 



A more detailed calculation for the MD cases would 
replace the ^-functions by Gaussians and the step func- 
tions by the corresponding error functions fl6| , with no 
essential change from the approximate results given here. 
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